#include "AnalysisUtils.hxx"
#include <TVector3.h>
#include <stdio.h>
#include <math.h>
#include <iostream>
#include <typeinfo>
#include "BaseDataClasses.hxx"
#include "DetectorDefinition.hxx"
#include <BeamBunching.hxx>
#include "PionInteractionSystematic.hxx"
#include "SystematicBoxUtils.hxx"

BeamBunching bunching;

//**************************************************
Int_t anaUtils::GetRunPeriod(Int_t run){
//**************************************************

    if      (             run<6000)   return 0;//run1  water                
    else if (run>=6462  && run<7664)  return 1;//run2  water                
    else if (run>=7665  && run<7754)  return 2;//run2  air                  
    else if (run>=8000  && run<8550)  return 3;//run3b air => mc associated to be the one of run2                  
    else if (run>=8550  && run<8800)  return 4;//run3c air                  
    else if (run>=8983  && run<9414)  return 5;//run4  water                
    else if (run>=9426  && run<9800)  return 6;//run4  air                  
    else if (run>=10252 && run<10334) return 7;//run5  ??               
    else if (run>=10334 && run<10519) return 8;//run5  water (ANTINU)       


    else if (run>=90110000 && run<=90119999) return 0;//run1 water          
    else if (run>=90210000 && run<=90219999) return 1;//run2 water          
    else if (run>=90200000 && run<=90209999) return 2;//run2 air            
    else if (run>=90300000 && run<=90300015) return 3;//run3 air separated based on the pot of data (run3b/(run3b+run3c)=0.13542             
    else if (run>=90300016 && run<=90300110) return 4;//run3 air separated based on the pot of data (run3b/(run3b+run3c)=0.13542            
    else if (run>=90410000 && run<=90419999) return 5;//run4 water          
    else if (run>=90400000 && run<=90409999) return 6;//run4 air            
    else if (run>=80510000 && run<=80519999) return 8;//run5 antinu-water

    else if (run>=80307000 && run<=80307999) return 8;//sand muons antinu
    else if (run>=90307000 && run<=90307999) return 4;//sand muons neut
    else if (run>=92307000 && run<=92307999) return 4;//sand muons old neut


    // genie
    else if (run>=91110000 && run<=91119999) return 0;//run1 water   
    else if (run>=91210000 && run<=91219999) return 1;//run2 water          
    else if (run>=91200000 && run<=91209999) return 2;//run2 air            
    else if (run>=91300000 && run<=91300015) return 3;//run3 air separated based on the pot of data (run3b/(run3b+run3c)=0.13542             
    else if (run>=91300016 && run<=91300110) return 4;//run3 air separated based on the pot of data (run3b/(run3b+run3c)=0.13542            
    else if (run>=91410000 && run<=91419999) return 5;//run4 water          
    else if (run>=91400000 && run<=91409999) return 6;//run4 air            
    else if (run>=81510000 && run<=81519999) return 8;//run5 antinu-water
    else return -1;

}


//**************************************************
int anaUtils::GetSandMode(int run){
//**************************************************

  if      (run>=80307000 && run<=80307999) return 0;	// antineutrino flux
  else if (run>=90307000 && run<=90307999) return 1;	// neutrino flux
  else if (run>=92307000 && run<=92307999) return 2;	// neutrino flux, old NEUT

  else return -1;
}

//*****************************************************************************
int anaUtils::GetFGDMichelElectrons(const AnaEventB& event, const SubDetId::SubDetEnum det, AnaFgdTimeBinB** arr, bool prod5Cut){
//*****************************************************************************  

    int nSel = 0;

    Int_t ibunch = event.Bunch; 
    Float_t selTime = bunching.GetBunchCentralTime(event, ibunch); 

    Float_t numWidthsReqInSorting = 4; 
    Float_t binBunchWidth = 25.; 
    Float_t binTimeOffsetToBunch = -13.;
    Float_t minRecBinQ = 0.; 
    Int_t minHitsFGD1 = 6;
    Int_t minHitsFGD2 = 5;

    if(prod5Cut){
        minRecBinQ = 200.0;
        minHitsFGD1 = 0;
        minHitsFGD2 = 0;
    }

    int index=-1;

    selTime += numWidthsReqInSorting*binBunchWidth;

    for( Int_t i = 0; i < event.nFgdTimeBins; i++ ){

        AnaFgdTimeBinB *FgdTimeBin = event.FgdTimeBins[i];

        if( det == SubDetId::kFGD1 && FgdTimeBin->NHits[0] == 0 ) continue; 
        if( det == SubDetId::kFGD2 && FgdTimeBin->NHits[1] == 0 ) continue;
        if( det == SubDetId::kFGD  && FgdTimeBin->NHits[0] == 0 && FgdTimeBin->NHits[1] == 0 ) continue;

        Float_t binTime = FgdTimeBin->MinTime - binTimeOffsetToBunch;

        bool inBunch = false;

        for( UInt_t k = 0 ; k < bunching.GetNBunches(event); k++ ){
            Float_t time =  bunching.GetBunchCentralTime(event, k); 
            if( TMath::Abs( binTime - time ) < numWidthsReqInSorting*binBunchWidth )  inBunch = true;
        }

        if( inBunch )  continue;

        //find bin with min raw charge and greater than the required number of hits
        if( binTime >= selTime && FgdTimeBin->RawChargeSum[0] > minRecBinQ && FgdTimeBin->NHits[0] > minHitsFGD1 &&  (det == SubDetId::kFGD || det == SubDetId::kFGD1) ){
            index = i;
        }
        if( binTime >= selTime && FgdTimeBin->RawChargeSum[1] > minRecBinQ && FgdTimeBin->NHits[1] > minHitsFGD2 && (det == SubDetId::kFGD || det == SubDetId::kFGD2) ){
            index = i;
        }

    }
    if(index > -1){
        arr[nSel] = event.FgdTimeBins[index];
        nSel++;
    }

    return nSel; 
}

//********************************************************************
std::vector<Float_t> anaUtils::ComputeWeights(const AnaEventB& event){
//********************************************************************

    std::vector<Float_t> weights;

    //  Only beam flux reweight for the moment
    weights.push_back(event.Weight);

    return weights;
}

//********************************************************************
void anaUtils::ClearArray(AnaTrackB** objArr, int nObj){
//********************************************************************

    for(int i = 0; i < nObj; ++i){
        objArr[i] = NULL;
    }
}

//********************************************************************
void anaUtils::ClearArray(AnaTrueTrackB** objArr, int nObj){
//********************************************************************

    for(int i = 0; i < nObj; ++i){
        objArr[i] = NULL;
    }
}

//********************************************************************
void anaUtils::FillArray(AnaTrackB** tgtArr, AnaTrackB** srcArr, int nObj){
//********************************************************************

    for(int i = 0; i < nObj; ++i){
        tgtArr[i] = srcArr[i];
    }
}

//********************************************************************
void anaUtils::FillArray(AnaTrueTrackB** tgtArr, AnaTrueTrackB** srcArr, int nObj){
//********************************************************************

    for(int i = 0; i < nObj; ++i){
        tgtArr[i] = srcArr[i];
    }
}

//********************************************************************
void anaUtils::InitializeEventSystematics(const AnaEventB& event, bool initialisePionSystematics,bool dofgd1,bool dofgd2){
//********************************************************************
  SystematicBoxUtils boxUtils(event,dofgd1,dofgd2);
  boxUtils.InitializeEventSystematics(event,initialisePionSystematics);
}
//********************************************************************
std::string anaUtils::GetSoftwareVersionFromPath(const std::string& path){
//********************************************************************

    std::string version="";
    size_t lastSlash = path.find_last_of("/");
    version = path.substr(lastSlash+1);

    return version;
}
//*****************************************************************************
bool anaUtils::IsFGD1Analysis(const AnaEventB& event){
//*****************************************************************************

  if(!event.Summary)              return false;
  if(!event.Summary->EventSample) return false;

  if((event.Summary->EventSample >= nd280Samples::kFGD1NuMuCC
      && event.Summary->EventSample<=nd280Samples::kFGD1NuMuCCOther) ||
     (event.Summary->EventSample >= nd280Samples::kFGD1AntiNuMuCC
      && event.Summary->EventSample<=nd280Samples::kFGD1AntiNuMuCCnQE)||
     (event.Summary->EventSample >= nd280Samples::kFGD1NuMuBkgInAntiNuModeCC
      && event.Summary->EventSample<=nd280Samples::kFGD1NuMuBkgInAntiNuModeCCnQE))
    return true;
  return false;
}
//*****************************************************************************
bool anaUtils::IsFGD2Analysis(const AnaEventB& event){
//*****************************************************************************
  if(!event.Summary)              return false;
  if(!event.Summary->EventSample) return false;

  if((event.Summary->EventSample >= nd280Samples::kFGD2NuMuCC
      && event.Summary->EventSample<=nd280Samples::kFGD2NuMuCCOther) || 
     (event.Summary->EventSample >= nd280Samples::kFGD2AntiNuMuCC 
      && event.Summary->EventSample<=nd280Samples::kFGD2AntiNuMuCCNTracks) )
    return true;
  return false;
}
